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We study the enzymatic degradation of an elastic fiber under tension using an anisotropic random- 
walk model, coupled with binding-unbinding reactions that weaken the fiber. The fiber is represented 
by a chain of elastic springs in series, surrounded by two layers of sites along which enzyme molecules 
can diffuse. Through numerical simulations we show that the fiber stiffness decreases exponentially 
with two distinct regimes. The time constant associated with the first regime decreases with in- 
creasing applied force, which is in agreement with published experimental data. In addition, a 
simple mean field calculation allows us to partition the time constant into geometrical, chemical 
and externally controllable factors, which is corroborated by the simulations. 
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The extracellular matrix (ECM), the biological struc- 
ture that supports cells, is composed of elastic fibers such 
as elastin and collagen. The complex organization of 
these fibers undergoes a continuous maintenance that re- 
uires the catalytic action of enzymes, called proteases 
1]. In diseases, such as pulmonary emphysema, tissue 
destruction is thought to be a consequence of the imbal- 
ance between protease and antiprotease activities leading 
to degradation of elastin fibers Q. Biological tissues in 
vivo are also under tension which may interfere with the 
enzymatic activity. Indeed, recent experiments show that 
mechanical stretch accelerates the rate of degradation of 
engineered ECM during elastase- induced digestion [3[. 

The elasticity of a single fiber depends on how its 
molecular constituents are organized. During digestion, 
the molecules in the fiber as well as the cross-links can 
be cleaved by enzymes causing fiber stiffness to decreases. 
Furthermore, following cleavage, an enzyme can unbind, 
diffuse, bind at a different location and cleave another 
molecule. This leads to the question: How are the dif- 
fusion and binding of the enzyme and the subsequent 
degradation of the fiber affected by the presence of an 
external mechanical force on the fiber? 

In this letter, we study the decay of stiffness of a sin- 
gle fiber under tension during enzymatic digestion using 
an anisotropic random walk model coupled with binding- 
unbinding reactions. To our knowledge this is the first 
investigation of the mechanical properties of a single fiber 
that takes into account the simultaneous effects of en- 
zyme diffusion, binding, cleaving and mechanical force. 

Several different diffusion-reaction models have been 
used to describe processes at the level of ECM, cell mem- 
branes, macromolecules and DNA 0, IE IE IS 0- In 
addition, the random walk is often used as a diffusion 
model that takes into account the morphological details 
of the system. Spring network models also provide a 
useful framework for analyzing the ch ang es in the me- 
chanical properties of ECM sheets [lfj, |ll| . Here we use 
a modified random walk to mimic anisotropic diffusion 
of enzyme particles along a fiber and study the digestion 



of the fiber represented by a chain of springs. 

Our model consists of a one-dimensional chain of N s 
linearly elastic springs in series representing an elastic 
fiber as in Fig. [TJi. The fiber is surrounded by two layers 
of sites along which particles representing enzymes can 
diffuse. Periodic boundary conditions are applied in the 
x direction. Both ends of the chain are subjected to a 
constant force F during digestion which mimics tension 
in the fiber. 

In order to simulate enzyme activity on the fiber, we 
begin with a chain having identical initial spring con- 
stants k{t = 0) = fco. The diffusion of enzyme is initi- 
ated by releasing a set of particles at random positions in 
the two layers. Each particle moves according to a set of 
probabilistic rules, controlling the diffusion and reaction 
processes: (i) pd is the probability for a particle to move 
right or left, parallel to the chain. This probability is 
associated with diffusion, (ii) p on is the probability for a 
particle to move up from the bottom or down from the 
top layer. This step represents an enzyme molecule bind- 
ing to a binding site on the fiber. Only one particle can 
be bound to a single spring at any time, (iii) p ff is the 
probability for a bound particle to move up or down to 
the top or bottom layer, respectively. This step is related 
to an enzyme molecule unbinding from the fiber. Once a 
bound particle unbinds, the local spring contant k is re- 
duced by a constant factor j, k — ► jk. The probabilities 
p on and pd are related by the constraint p on + 2pd = 1. 

Local anisotropy is introduced through the probability 
p on which depends on the local spring constant k, 



(1) 



where Ap is the initial anisotropy and A is a character- 
istic length. Eq. (fT]) expresses the fact that when the 
enzyme cleaves the fiber, the local k decreases by a small 
amount. Since F is constant, the local stretch increases 
and either more binding sites appear or the binding of 
an enzyme becomes easier which we represent by an in- 
crease in the local p on . This introduces anisotropy in the 
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FIG. 1: (a) Schematic diagram of the chain of springs and 
sites used in our model. The binding sites on the springs and 
the two layers of sites are represented by small and big open 
circles, respectively. The set of enzyme particles are shown as 
filled black circles. The particle at the botton layer can move 
up, left, or right while the particle at the top can move down, 
left, or right. The particle on the spring can move only up or 
down, (b) The binding probability p on defined by Eq. ([TJ and 
the diffusion probability pd as a function of the number of 
visits n at a fixed site. Lines of different styles correspond to 
Ap = 0.10 (solid lines), Ap = 0.16 (dotted lines), Ap = 0.20 
(dashed lines) and Ap — 0.24 (dash-dotted lines). The lines 
above and below the horizontal dashed line at Ap = 0.333 
corespond to pd and p on , respectively. The vertical dashed 
line represents the region where the isotropic behavior with 
p — 0.333 is reached. 



particle movement and makes the enzyme activity de- 
pendent on the local k. The parameter Ap is related to 
difficulty of an enzyme molecule to reach a binding site 
which depends on the surface roughness of the fiber. In 
Fig. QJ we can see that at the beginning of the diffusion, 
the number of times a spring has been visited by parti- 
cles is ti w 1 and the initial values of both p on and p,i 
depend on Ap. As the diffusion progresses, p on increases 
slowly. Around n = 100, p on increases significantly until 
it approaches the isotropic value case p on = pd = 1/3, 
where the rate of increase is reduced. When n 800, 
the diffusion reaches the isotropic regime in which both 
diffusion and binding are equally likely. We assume that 
in this regime, locally the fiber is at its unfolding limit in 
that the number of binding sites remains constant. 

The probability p Q ff is determined by the molecular 
properties of the specific enzyme and its substrate and is 
the same for each spring, for all times. A particle may re- 
main bound for more than one time step, with probability 
1 — PoB- We assume that cleavage occurs during unbind- 
ing so that k is reduced only when the particle unbinds 
from the spring. Thus, k for each spring is a function 



of time because it decreases as the spring is repeatedly 
visited by particles. Since k decreases by the constant 
factor 7, we can write the relation k(t) = 7™fco, where n 
is the total number of visits by time t. The stiffness K (t) 
of the fiber is calculated as the equivalent stiffness of all 
N s springs connected in series, K(t) — l/J^it^C*)] -1, 
Next, we study the evolution of K(t) for differ- 
ent sets of parameters. We use a chain composed of 
N s = 10 4 springs, different numbers of particles N p = 
256, 512, 1024 and different values of the external force F 
within the interval [0.1,2.5]. The N p and F are related 
to the experimentally controllable macroscopic parame- 
ters of the digestion process. At the microscopic level, we 
vary p Q g between [0.1, 1.0] while Ap is chosen from the in- 
terval [0.10,0.24]. Additionally, there are three constant 
parameters: A = 0.10, 7 = 0.995 and ko = 1. We obtain 
the time course of K(t) for i = 2x 10 5 time steps, where 
at each time step, we attempt to move all N p particles 
in the system. We repeat the digestion simulation 500 
times with different realizations and average K{i) over 
all runs. 

The results for (K(t)) are plotted on a log- linear scale 
in Fig. [21 In all cases K(t) shows two distinct exponen- 
tially decreasing regimes with time constants T± and Ti 
separated by a crossover region around t x . For fixed N p , 
T\ decreases monotonically as F increases. The t x also 
decreases as a function of F leading to a faster overall 
degradation of the fiber. The situation is similar when 
N p is increased. 

To characterize the changes in the microscopic prop- 
erties of the fiber, we calculate the standard deviation 
((jfe) of all spring constants at a fixed time and average 
them over all runs in Fig. [3jz. Initially, for t <C t±, (fife) 
increases quickly, which is not influenced by F. When 
t « 2000, F starts to affect the binding process accord- 




FIG. 2: Log-linear plot of the averaged stiffness (K) as a func- 
tion of diffusion time t for different values of F = 0.1, 0.5, 2.5 
and N p = 256, 512 with Ap = 0.20 and p off = 0.5. The black 
solid line segments at the beginning and end of the simula- 
tions represent exponential fits to estimate the value of the 
time constants Ti and T%, respectively. 
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FIG. 3: (a) The standard deviation (<jfc) of local spring 
constants k as a function of time for three values of F = 
0.1,0.5,2.5 and 7V P = 512. In (&), we plot the distribution 
of spring constants P(k) at time points ti (i = 1,2,3) as 
indicated in (a) for F = 0.1. In both graphs, we use the 
parameters Ap =0.20 and p g = 0.5. 
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FIG. 4: Comparison of the mean field calculations and the 
numerical simulations. Time constants Ti and Ti are plotted 
as functions of F and p « for N p = 256 (circles), N p = 512 
(squares) and N p = 1024 (triangles). The symbols corespond 
to the numerical simulations and the solid lines are obtained 
from Eq. (J5J). In (a) and (b), Ti and T2 as a function of F for 
a constant p a g = 0.5. In (c) and (d), Ti and Tb as a function 
of Pofi for F = 1.0. In all graphs we used Ap — 0.20. 



where we also assume that k m ~ (fc). Thus, from Eqs. @ 
and ([3]), the change AK = K(ti + i) — K(ti) in the total 
stiffness is written as 



ing to Eq. ([I]) and hence (<t/c) increases faster for higher 
F. When t is around t 2l (ffc) reaches its maximun fol- 
lowed by a slow decay. The behavior can be confirmed 
by looking at the spring constant distributions P(k) in 
Fig. [3f). The width of P(k) has a maximum at t 2 . For 
t < ti and t > P(k) becomes narrower. Also, the peak 
of P(k) decreases with increasing time. Similar behavior 
is observed for higher iV p ; however, the maximum (<7fc) 
decreases when N p increases (data not shown). 

To gain insight into the exponential decay of fiber stiff- 
ness, we carry out a simple mean field calculation. The 
local k at time t depends on the number of times the 
spring has been visited. We define U as the time cor- 
responding to the i-th unbinding event along the entire 
chain. Thus, K at time ti, can be written as 



^?=1 k. 



(fc) 
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assuming (1/fc) ~ Notice that K remains constant 

for ti ^ t < At time ii+i, an unbinding event occurs 
at spring m and the corresponding k m is reduced to 7fc m . 
The new value of K(ti + i) is 



AK = K(t i+1 ) - K(ti 
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Next, we consider the average waiting time (t) between 
two unbinding events. During digestion, the number of 
particles nn that remains bound on the fiber changes, 
but ns is related to the number of free particles np since 
ris + np = N p . The rate of change of ns is the dif- 
ference between the average binding rate p on np and 
the average unbinding rate p Q g n^ which can be ex- 
pressed as dns/dt = p on np — p g ns- Assuming that 
nB has reached a steady state dns/dt = 0, we obtain 

n B =PonNp/ (p on +p s)- 

Since the unbinding probability is per unit time, (t) 
can now be expressed as 
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Finally, we can establish a link between the two pro- 
cesses involved in enzymatic digestion. From Eq. K 
is reduced by AK during the interval At — — ti. 
We thus approximate the derivative of K by the discrete 
change AK during the interval At ~ (t) as 
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This equation can be solved assuming that p on is ap- 
proximately constant during one time step. The result is 
given by Kit) — e~*/ T , where 
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The above expression represents different aspects of the 
digestion process. The first term which involves 7, is 
related to the geometry of the fiber, the average num- 
ber of molecules in parallel. The second term is related 
to the specific enzyme activity at the microscopic level, 
the binding and unbinding processes. The third term 
is essentially the enzyme concentration that is an exter- 
nal control parameter. Note that the time constant T 
is symmetric in the binding and unbinding probabilities. 
Thus, if we assume it is the unbinding process that de- 
pends on the external force, the results for the stiffness 
degradation will be identical. 

We next analyze the asymptotic limits of Eq. ([7]). We 
summarize the results as follows 
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To compare the results of the numerical simulations 
with the analytical calculation in the asymptotic limits, 
we analyze the time course of stiffness by calculating the 
time constants in Fig. [2] for the two different regimes of 
exponential behavior. We obtain exponential fits in non- 
overlapping windows with fixed size At = 5000 along 
the stiffness curve. The T\ and T2 are then extracted 
as a function of F and p Q s in the regions Kit) ~ 1 and 
Kit) <C 1, respectively. The results are compared to 
those obtained from Eq. {5} in Fig. 0] Generally, the re- 
sults confirm the agreement of T\ and Ti measured in the 
simulations and in the analytical calculations as a func- 
tion of F and p g. The T\ decreases with increasing F 
while it diverges when p Q g becomes less than 1/3. Notice 
also that for fixed values of F and p s, T\ decreases with 
increasing enzyme concentration. Interestingly, F has lit- 
tle effect on T^. This can also be seen in Fig. [5] which 
shows that after the crossover region all curves follow the 



same exponential decay for different F at a fixed value of 
N p , but decays slower as N p decreases. The difference in 
T2 between analytical and numerical calculations for low 
N p is because K, in this limit, takes significantly more 
time to reaches the second exponential regime. To con- 
firm this, we run an additional simulation for N p = 256 
for 3 x 10 5 time steps and find that T2 increases by about 
15%, approaching the analytical results. Also for a fixed 
F and p a s, the T\ increases as Ap increases and T 2 re- 
mains approximately constant as a function of Ap (data 
not shown). 

In summary, we have presented a model for the en- 
zymatic digestion of an elastic fiber under tension. We 
have shown that the total stiffness decreases exponen- 
tially with two different regimes separated by a crossover 
region. While the first exponential regime has been found 
experimentally [|[, to our knowledge, the second one 
has not been measured and remains a prediction of the 
model. Each regime can be associated with an average 
value of the stiffness along the fiber during the particle 
diffusion process. In the first regime, the stiffness is domi- 
nated by the average local initial stiffness values whereas 
in the second regime, the local stiffness has decreased 
significantly and almost uniformly throughout the fiber. 
In the crossover region, the stiffness is controlled by a 
wide distribution of local stiffness values. The time con- 
stant T\ displays a strong dependence on both F and N p . 
This result is in agreement with the experimental results 
3 where it was reported they show that application of 
static mechanical forces accelerates the digestion-induced 
breakdown of ECM sheets. Also, we have performed an- 
alytical calculations that confirm the presence of two dif- 
ferent regimes. These calculations also show how the 
time constant can be partitioned into geometrical, chem- 
ical and externally controllable factors. Furthermore, in 
the first regime, we expect that the fiber does not reach 
the failure limit but after the crossover the decrease in 
the stiffness is more intense and failure is likely to occur. 
These results can help better understand natural growth 
and maintenance of the ECM as well as diseases in which 
enzyme concentrations and/or mechanical forces become 
abnormal. 
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